A mathematically rigorous algorithm to define, compute and assess relevance of the probable dissociation constants in characterizing a biochemical network

Metabolism results from enzymatic- and non-enzymatic interactions of several molecules, is easily parameterized with the dissociation constant and occurs via biochemical networks. The dissociation constant is an empirically determined parameter and cannot be used directly to investigate in silico models of biochemical networks. Here, we develop and present an algorithm to define, compute and assess the relevance of the probable dissociation constant for every reaction of a biochemical network. The reactants and reactions of this network are modelled by a stoichiometry number matrix. The algorithm computes the null space and then serially generates subspaces by combinatorially summing the spanning vectors that are non-trivial and unique. This is done until the terms of each row either monotonically diverge or form an alternating sequence whose terms can be partitioned into subsets with almost the same number of oppositely signed terms. For a selected null space-generated subspace the algorithm utilizes several statistical and mathematical descriptors to select and bin terms from each row into distinct outcome-specific subsets. The terms of each subset are summed, mapped to the real-valued open interval \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\left(0,\infty \right)$$\end{document}0,∞ and used to populate a reaction-specific outcome vector. The p1-norm for this vector is then the probable dissociation constant for this reaction. These steps are continued until every reaction of a modelled network is unambiguously annotated. The assertions presented are complemented by computational studies of a biochemical network for aerobic glycolysis. The fundamental premise of this work is that every row of a null space-generated subspace is a valid reaction and can therefore, be modelled as a reaction-specific sequence vector with a dimension that corresponds to the cardinality of the subspace after excluding all trivial- and redundant-vectors. A major finding of this study is that the row-wise sum or the sum of the terms contained in each reaction-specific sequence vector is mapped unambiguously to a positive real number. This means that the probable dissociation constants, for all reactions, can be directly computed from the stoichiometry number matrix and are suitable indicators of outcome for every reaction of the modelled biochemical network. Additionally, we find that the unambiguous annotation for a biochemical network will require a minimum number of iterations and will determine computational complexity.


Terminology and concepts pertinent to reaction kinetics for use by a stoichiometry number-based model of a biochemical network
The formation of a complex is common to enzyme-mediated catalysis as well as non-catalytic cycles of associations and disassociations.We will utilize the term "reaction" throughout the manuscript to indicate this, We extend this schema to the encompassing pathway and utilize the phrase "biochemical network" to represent a set of reactants/products that participate in a set of reactions (Def.(2)).We utilize the stoichiometry number for each reactant/product to indicate the presence of a change or lack thereof, that will be brought about by a reaction (Def.(3)).Consider the reaction (r) with j-indexed J-reactants and with rate (R r ) .We will define the rate constant ( r ) as, The intracellular environment is extremely complex and will support several forms of reactions between the reactants/products whilst precluding others 28,[33][34][35] .For example, the intracellular milieu or cytosol is regarded as non-Newtonian and imposes constraints on the type of reactions that may occur along with changes that effect reaction rate 34,35 .These alterations in local viscosity can occur in response to acute stressors such as fluctuations in temperature and hydrogen ion concentration.The modifications have a biochemical (trehalose, glycogen) or biophysical (cytoskeleton) basis and may function to modulate diffusion and diffusion-mediated chemical reactions in vivo or in silico [34][35][36][37] .The physiological manifestations of these will include reduced motion of macromolecules, increased solubility and phagocyte movement towards a chemotactic and noxious stimulus [36][37][38] .In order to retain this complexity we will use a reaction vector to model the contribution of each reactant/product.For the set of contributing j-indexed J-reactants or products, Definition 5 (D5) The rate constant for a reaction ( r ) is a strictly positive real scalar, Although a reaction may have several outcomes, we will only consider the forward-, reverse-and equivalentoutcomes in this study.Here, we define reaction "equivalence" as the uncertainty with regards to the preferred outcome for a reaction (Def.(6)).This may be due to the paucity of data on the reactants/products and/or the prevailing biophysical state of the intracellular environment where the reaction is expected to occur.In contrast, the reversibility of a reaction (forward, reverse) is a well characterized outcome.

The dissociation constant as the parameter of choice to investigate a biochemical network
Parameterizing a chemical reaction is non-trivial and is done empirically, derived theoretically or inferred from data 23,29,33,39 .For a non-enzymatic reaction the dissociation (Kd) -and association (Ka)-constants are used to describe reversibility (forward, reverse) at or near equilibrium 4,5,33,40 ,
The dissociation constant also allows us to describe a "non-productive" reaction which could occur when two or more reactants/products associate, form a complex and remain associated (Def.(7)).Unlike non-reacting reactants/products, the reaction commences but does not terminate.This can happen when the product of an enzyme forms a covalent linkage with the active site residues and inactivates it irreversibly [41][42][43] .A non-productive reaction is also observed when an inhibitor binds preferentially to the bound form of an enzyme, i.e., enzyme-substrate complex, rather than the free enzyme.This unique occurrence is referred to as "uncompetitive"-inhibition and is routinely observed in the tumor protein (TP)53-mediated inhibition of Glucose-6-phosphate dehydrogenase (G6PD; EC1.1.1.49) as a well-supported mechanism for tumor suppression 43 .For a non-enzymatic reaction this implies an exceptionally high binding affinity between the reactants and/or products with the consequence that these become untenable as independent reactants/products [41][42][43] .
The dissociation constant is a non-negative real number (Kd ∈ R ∩ (0, ∞) ) , can be empirically determined (surface plasma resonance, fluorescence, spectroscopic methods) or theoretically derived, and is a reliable index of reaction outcome 33,[44][45][46][47] , Since the dissociation constant is a ratio, we can use it as a functional metric to compare the rates at which each reaction of a biochemical network takes place as well as the rates of the forward-and reverse-components of a reversible reaction, The dissociation constant is versatile, informative and biochemically relevant and is therefore, the parameter of choice to investigate a biochemical network in our study.Here, we will derive a numerical measure directly from the stoichiometry number matrix of a biochemical network which we define as the "probable dissociation constant" (Def.( 8)).

Stoichiometry number-based model of a biochemical network
We now model a biochemical network p as the sparse stoichiometry number matrix S p of J × I dimensions S p ⊂ Z J×I .Briefly, each column of this matrix is an i-indexed i = 1, 2 . . .I collection of non-trivial reaction vectors r 1 r 2 . . .r I ⊂ Z J .Each row of this matrix is contributed to by j-indexed j = 1, 2 . . .J m-stoichiometry numbers m ji ∈ Z of J-reactants/products across all I, Consider the forward reactions: Consider the reverse reactions: Consider the equivalent reactions: Vol.:(0123456789) www.nature.com/scientificreports/ The null space-generated subspace can be used to compute the probable dissociation constant for a reaction For a biochemical network p under the assumption of equilibrium, a subspace of the null space of a stoichi- ometry number matrix V ⊆ N S p is defined by the null space spanning vectors and those that result from their combinatorial sums (Fig. 1), Clearly, the cardinality of a null space-generated subspace will increase with the number of cycles of combinatorial summations which we define as u-iterations (Def.(12)).The comprehensive subspace for the uth-iteration (V u ) is then defined by combinatorially summing the vectors from the (u − 1)th-iteration (V u−1 ) .Using this notation we can easily see that, Definition 13 A null space-generated subspace for the uth-iteration is comprehensive and defined as, In order to partition H u into subsets we will utilize the following definitions and notation, where, 2) where, N S p = Null space spanning vectors for a given stoichiometry number matrix Figure 1. Outline of the algorithm to compute the probable dissociation constant and assign outcome to every reaction of a biochemical network.A biochemical network is modelled in terms of the stoichiometry numbers for a finite number of reactants/products and the accompanying reactions.A null space-generated subspace comprises only unique and non-trivial vectors which are iteratively, recursively and combinatorially summed.The reaction-specific sequence vector that is formed after a finite number of iterations will comprise terms that do not converge and are easily partitioned into distinct subsets.Statistical and mathematical descriptors for these terms (mean, standard deviation, greatest lower bound, least upper bound) are used to select and bin terms for all three predicted outcomes of a reaction.The formed reaction-specific outcome vector has a p1-norm which is the probable dissociation constant for a reaction.This process is repeated until every reaction of the modelled biochemical network can be unambiguously assigned.I , Number of reactions for a modelled biochemical network; J , Total number of reactants of modelled biochemical network; S p , reaction-centric and user-defined stoichiometry number matrix for a biochemical network; A u ; u-iteration specific subspace of the null space of the stoichiometry number matrix for the modelled biochemical network; {F i , B i , E i } , Subsets to bin kth-term of an ith-reaction-specific sequence vector and (u = u > M)th-iteration; g(.) ; Map to assign the sum of F i -, B i -and E i -subsets to a strictly positive real number of an ith-reaction-specific sequence vector and (u = u > M)th-iteration; T, C, D; theorem, corollary, definition; Vol.:(0123456789) Interestingly, we can use T1 to compute the lower bounds for a biochemical network in terms of the reactants/ products (J) and reaction vectors I .Let us also assume that a single master reaction is an improbable event for a functional biochemical network given the complexity of the intracellular milieu, We see that in order to remain biochemically relevant the modelled network must have at least, The corresponding number of reaction vectors will be,

Computational tools to calculate the probable dissociation constant for every reaction of a biochemical network
The biochemical relevance and suitability of the probable dissociation constants for a biochemical network is illustrated by characterizing (constructed, presented, analysed) a biochemical network for AG.The R-package "ReDirection" will be used to assess the contribution of reduced-and oxidized-forms of Nicotinamide adenosine dinucleotide phosphate in regulating the activity of the Pyruvate dehydrogenase complex (PDC) 48 ."ReDirection" is freely available and can be accessed without a login id from the comprehensive R archive network (CRAN; https:// cran.r-project.org/ packa ge= ReDir ection) 48 .
Briefly, "ReDirection" will check and modify, if necessary, the stoichiometry number matrix for a biochemical network 48 .The matrix, along with a logical argument (TRUE, FALSE) which indicates whether the reactions are to be regarded as rows or columns, are entered as arguments for the function "calculate_reaction_vector" 48 ."ReDirection" will exclude all linear dependent vectors and then recheck the modified matrix for compliance with network-specific indices such as the minimum number of reactants/products, minimum number of reactions and a numerical difference of at least two in favour of the number of reactions 48 .These are necessary arguments and "ReDirection" will not process the stoichiometry number matrix if these are violated.In case there are no deficiencies, "ReDirection" will compute the null space and generate subspaces serially by combinatorially summing all non-redundant and non-trivial vectors 48 ."ReDirection" will continue these steps for a finite number of iterations until a subspace is found where the lower bounds for each sequence of row terms are unambiguous and well-defined."ReDirection" will then check and bin every selected term to output-specific subsets (forward, reverse, equivalent), compute the linear sum for each subset, map these to strictly positive real numbers and thence populate a reaction-specific output vector 48 .The p1-norm of the latter is denoted as the probable dissociation constant for a reaction 48 ."ReDirection" will repeat this until every reaction is annotated unambiguously 48 .The output will comprise a list of reactions with their computed probable dissociation constants and the predicted outcome (forward, reverse, equivalent) 48

Computational studies to highlight and assess biological relevance of the probable dissociation constants
Although the probable dissociation constants that "ReDirection" computes is specific for a biochemical network, there are several generic metrics which can be derived such as the proportion and distribution of matched reactions (equivalent, non-equivalent).These can be used to indicate the flux of reactants/products in a particular direction within a network and compare the modelled biochemical network under variable intracellular conditions.
Aerobic glycolysis occurs when Glucose is converted to Lactate in the presence of non-limiting molecular dioxygen 49,50 .Although purported initially as a plausible mechanism to explain the atypical metabolism of tumorous cells, the phenomenon is routinely observed in rapidly proliferating non-tumorous (enterocytes, hematopoietic stem cells), quiescent (fibroblasts) and endothelial cells of newly forming vasculature 50 .AG offers plausible explanations for mechanisms for innate immune cell memory, angiogenesis and macrophage polarization [50][51][52][53] .For example, a common observation in pro-inflammatory polarized macrophages is a break in Kreb's cycle at the level of citrate [54][55][56][57] .However, there is no clarity on the genesis and reversal of this outcome in the presence an intact and functioning Kreb's cycle.Conversely, the roles of the bifunctional enzyme phosphofructokinase-2/ fructose-2, 6-bisphosphatase 3 (PFKFBP3) in the flux of pyruvate prior to its entry into the mitochondria have been implicated in several known AG-pathways 58,59 .Central to these hypotheses is a role for the PDC and its Vol:.( 1234567890)

Results and discussion
The steady-state for a system is described as the absence of apparent change in the numbers of the reactants/ products.For a biochemical network this is usually observed when the rates of formation and utilization of the reactants/products that comprise it are assumed to be equal.These conditions, when they occur within the cell are likely to be biochemical in origin and will include the effects of molecular redundancy (alternate splicing, duplicated genes, pseudogenes), feedback (positive, negative), compartmentalization and shared reactants 22,31 .

Null space generated subspaces with reduced cardinality may improve time to complete annotation
Although the comprehensive null space-generated subspace for a modelled biochemical network is desirable, the vectors that comprise this may be trivial and redundant.This will lead to greater computational complexity and thence delay the time needed to completely annotate every reaction of a biochemical network 48 .It is therefore, imperative that these vectors are identified and excluded.The resulting subspace (A u ) is of reduced cardinality and is defined as, The following corollaries about the cardinality of A u hold and are easily established.

Corollary 1 (C1)
The expected cardinality of null space generated subspace for the matrix of stoichiometric numbers of the reactants/products of a biochemical network is always greater than the computed cardinality, Corollary 2 (C2; without proof) The cardinality of an arbitrary and recursively generated subspace of the null space of the matrix of stoichiometry numbers of the reactants/products for a biochemical network is always greater than the nullity of the null space,

Corollary 3 (C3; without proof) The cardinality of recursively generated subspaces of the null space of the matrix of stoichiometric numbers of the reactants/products of a biochemical network forms a monotonically increasing sequence,
Since the time taken to annotate every reaction of a biochemical network will depend on the combinatorial sums of all non-redundant and non-trivial null space generated subspace vectors, the cardinality of this subspace, at each iteration is an important determinant of computational complexity.

Corollary 4 (C4) The computational complexity per iteration in exponential-time (T u ) is,
These results suggest that as the number of null space generated subspace vectors increases there will be an increase in the number of iterations needed to unambiguously annotate a reaction.This in turn will result in an increased cardinality of each null space generated subspace.A corresponding increase in run-time will then be needed to completely annotate every reaction of a biochemical network 48 .
S p ⊂ Z J×I = Stoichiometry number matrix for the biochemical network N S p = Null space of S p , C1 and Def .( 13)] (42) Vol.:(0123456789) www.nature.com/scientificreports/Rationale to utilize the row-sum of a null space-generated subspace as the domain to compute the probable dissociation constant for a reaction The vectors that comprise each null space-generated subspace represent a biochemical network at equilibrium where each row represents an ith-reaction and each numerical value belongs to R ∩ (−∞, ∞) .However, since the dissociation constant is a non-negative real number we seek a map to the set of strictly positive real numbers such that in the absence of empirical data we can directly compute the probable dissociation constant from the null space of the stoichiometry number matrix for a biochemical network, Let us partition the real part of the numerical values of any non-trivial and non-redundant null space generated subspace vector, Combining these, This means we can only map the following values unambiguously, For v ik ≈ 0 however, this mapping schema for an "equivalent"-reaction will not suffice since, We can resolve this by treating each row of a null space generated subspace (A) as a finite series, If we do this iteratively (u) , recursively and across all columns of a null space-generated subspace (k = 1, 2 . . .K) and for all reactions (∀i) , we obtain the following expression, Theorem 2 (T2) The numerical values that comprise each row of a null space generated subspace (A u ) with posi- tive-and negative-terms can be rewritten as an alternating sequence, Theorem 3 (T3) The sum for the finite series formed by each row of a null space generated subspace can be mapped from R ∩ (−∞, ∞) to R ∩ (0, ∞) for the ith-reaction of the uth-iteration, Def.( 16) where, f := Index for a putative forward outcome predicting term Def.( 17) where, where, Vol:.( 1234567890) www.nature.com/scientificreports/ We can now utilize φ u i to unambiguously compute the probable dissociation constant and assign an outcome (forward, reverse, equivalent) to the i th -reaction and thence to every reaction of a biochemical network.

Algorithm to compute the probable dissociation constant for every reaction of a biochemical network
An R-package implementation of the basic ideas that underline the presented algorithm is available for a userdefined biochemical network 48 .However, the mathematical basis for these findings and/or observations have not been addressed.Here, we rigorously investigate some of the several fundamental steps in the computation, usage and utility of the probable dissociation constant for every reaction of a biochemical network.

Schema to screen combinatorially summed vectors for a null space-generated subspace
For each u-iteration the null space-generated subspace will be screened for trivial and redundant vectors and subsequently excluded.The non-redundant and non-trivial vectors that remain will constitute the cardinality for the u-iteration specific null space-generated subspace (Step 1, Fig. 1).The following assumptions will be valid for the presented schema, Rewriting the finite set of identical null space-generated subspace vectors for the u th -iteration H u in terms of a finite collection of non-overlapping subsets (�), Let us now define a representative vector from each subset as, Since there can be one and only one such vector for each subset, we can write, The finite series formed is then numerically equal to the number of subsets that H u is partitioned into, where, We will now reassign these vectors to and exclude these from Vol.:(0123456789) www.nature.com/scientificreports/These vectors will now be assigned to H which is the subset for all unique null space-generated subspace vectors.The incremented cardinality for the uth-iteration of this subset H u is, Reciprocally, H u will be correspondingly deficient in these vectors and can be written as, The reduced cardinality for the u th -iteration of H u is,

Characterizing reaction-specific sequence vectors from the rows of selected null space-generated subspaces
We have already seen that after a finite number of iterations of combinatorial summing and vector exclusion, null space-generated subspaces with reduced cardinalities will result (A u ).Since the terms from the ith-row when taken together, i.e., across all columns of the comprising vectors (cols(A u )) constitute the ith-reaction, we observe that #cols(A u )-terms can be partitioned into distinct subsets.We define u = M ∈ N formally as the lower bound for these iterations as (Step 3; Fig. 1).
The numerical values for M are easily computed from the limit of the fractional optimization of two decoupled real-valued functions when evaluated as a single ratio, Since the algorithm is based on combinatorial summations, the constituent real-valued terms of a reactionspecific sequence vector will be monotonic and their linear sum will be non-zero.This implies that convergence will occur very quickly (forward, reverse).If the linear sum of the terms of the reaction-specific sequence vector is zero, this will indicate that the sequence has an almost equal number of oppositely signed terms and will not converge (equivalent).However, here too, the use of fractional optimization, decoupling and a single ratio to minimize will result in convergence and the occurrence of a limit, where, Now that we have established the minimum number of iterations, we will proceed to define, populate and process the reaction-specific sequence vector for every reaction of the modelled biochemical network.We will use a simplified notation hereafter without loss of generality to indicate and annotate, We denote an ith-reaction-specific sequence vector sequence a i T with real values which for the uth-iteration whose sum is φ i in accordance with T2 (Step 4; Fig. 1), Once the ith-reaction specific sequence vector is defined we will characterize it using standard statistical descriptors (Step 5; Fig. 1), where where, u = 1, 2 . . .

Rule-based population of the outcome-specific subsets for every reaction of a biochemical network
We will first define the criteria to populate the F i -and B i -subsets from the reaction-specific sequence vector for the ith-reaction of a (u = U > M)th-iteration, Corollary 5 (C5) If the cardinality of the reaction-specific F i -and B i -subsets is finite we can identify a unique subset (E i ) with terms that form an alternating sequence after u = U > M-iterations,

Corollary 6 (C6, without proof)
The sum of the terms of the subset that corresponds to the "equivalent"-outcome is numerically equal to zero after a large u = U > M, where, where, and, e, f, g = 1, 2 . . .E, F, G ∈ N (87.We have already identified a generic mapping schema where the lower-and upper-bounds of the reaction-specific sequence vector that comprise the real-valued terms of the null space generated subspace can be mapped to the open interval R ∩ (0, ∞) (T2, T3; C5, C6) (Steps 6-9; Fig. 1).We will now define these maps formally for each outcome-specific subset, The p1-norm of the reaction-specific outcome-vector is the probable dissociation constant for that reaction Since every reaction can have these possible outcomes, we can define for the ith-reaction of the uth-iteration where u = U > M , a reaction-specific outcome vector, The p1-norm describes the final outcome of the ith-reaction after a finite number of u = U > M-iterations (Step 10; Fig. 1), where for a large u and, e, f, g = 1, 2 . . .
where, e, f, g = 1, 2 . . .E, F, G ∈ N i = 1, 2 . . .I Forward : a) For the non empty subset F i with cardinality #F i and sum φ i F , g : where, b) For the non empty subsets where, Reverse : For the subset B with cardinality #B i and sum φ iB , g : where, y iB = e φi B (98)

Equivalent :
For the non empty subset E with cardinality #E i and sum φ iE , g : where, Vol.:(0123456789) www.nature.com/scientificreports/ We can now summarize these assertions and unambiguously annotate every reaction of a biochemical network as, We can easily see that the probable dissociation constant (η i ) for the ith-reaction of the modelled biochemical network can be equated to the true dissociation constant (Kd i ), The algorithm is then iteratively applied to all other rows of a null space-generated subspace, i.e., reactions of the modelled biochemical network i = 1, 2 . . .I and the probable dissociation constants are computed for ∀i .If there is ambiguity in the annotation, the combinatorial summation and subsequent steps are recursively repeated 48 .Complete annotation for every reaction of a biochemical network will occur if and only if the annotation is unambiguous.A biochemical network is only regarded as being completely annotated if and only if every reaction that comprises the network is annotated unambiguously 48 .

Biochemical relevance and suitability of the probable dissociation constant to characterize a biochemical network of aerobic glycolysis
The biochemical network for AG I AG = 20; η 1 − η 20 comprises single ( I AG = 14 )-and multi ( I AG = 5)-step reactions in the cytosolic (c) and mitochondrial (m) compartments (Fig. 2) (Supplementary Text 2).Here, the multi-step reactions include the hexose monophosphate shunt (r 2 , r 4 ) , glycolysis (r 5 ) , Kreb's tricarboxylic acid (TCA; r 12 ), oxaloacetate-malate shuttle (r 16 ) and Cori's cycle (r 20 ) .The modelled biochemical network also includes transport reactions for Pyruvate (r 8 ) , Citrate (r 14 ) and mitochondrial Phosphoenolpyruvate (r 17 ) .The computed probable dissociation constants favour the export of oxaloacetate into the cytosol directly or through phosphophenolpyruvate ({η 10 , η 13 , η 16 , η 17 } > 40) thereby indicating a flux towards gluconeogenesis.The steps towards glucose 6-phosphate are supported by observing the probable dissociation constants for the following conversions: phosphoenolpyruvate → glyceraldehyde 3-phosphate (η 5 ≈ 0.017) and glyceraldehyde 3-phosphate → glucose 6-phosphate (η 3 ≈ 0.03) (Fig. 2).We can easily infer from the reverse reaction pyruvate → acetyly- CoA (η 9 ≈ 0.03) that the TCA-cycle is interrupted (Fig. 2).These findings suggest that contrary to available literature, the TCA-cycle is predisposed to being perturbed and its restoration or "intactness" is actively maintained for metabolic homeostasis or oxidative phosphorylation.These observations are perfectly logical given that the TCA-cycle is anaplerotic with several of its intermediates participating in collateral pathways.For AG to occur it is necessary that the TCA-cycle along with oxidative phosphorylation is interrupted in favour of an increased cell mass via the hexose monophosphate shunt pathway (DNA/RNA synthesis, NADPH) and fatty acid synthesis via acetyl-CoA (Fig. 2).
These data also offer several interesting insights into the metabolic crosstalk that may occur within cells when metabolizing Glucose to Lactate in the presence of adequate molecular dioxygen (Fig. 2).There is a significantly larger fraction 85%; I AG = 17 of non-equivalent reactions with forward 35%; I AG = 7 -and reverse 50%; I AG = 10 -reactions as compared to the equivalent 15%; I AG = 3 (Fig. 2).This suggests that the prevail- ing biochemical and physicochemical conditions where AG is expected to occur is, in fact an important determinant of whether AG occurs at all.Conversely, the paucity of equivalent reactions in the biochemical network for AG suggests lack of regulatory reactions.The ratio of NAD(P)H to NAD(P) + along with Acetyl-CoA to CoA and ATP to ADP is an important regulator of the Pyruvate dehydrogenase complex (PDC) and is likely a major determinant for the metabolic shift observed to implement reversible AG.

Limitations of the presented algorithm and relevance of fractional derivatives in modelling and analysing complex biochemical systems
The biomedical relevance in computing the probable dissociation constant notwithstanding, the presented algorithm is dependent on the stoichiometry or integer numbers to model change to the molecular species that comprise the modelled biochemical network.This effect is partially offset by computing the null space, although the assumption that the modelled biochemical network is at equilibrium is yet another presumption.A further limitation of the presented algorithm is that it is enumeration-based and will therefore, be intractable for large biochemical networks.Furthermore, since the computations are entirely dependent on the null space, the initial nullity is also likely to effect the real time needed to annotate every reaction of a modelled biochemical network.This will mean that for smaller networks too, the time to complete annotation may also not be possible 48 .
The probable dissociation constant is a numerical measure that is computed from the null space for a stoichiometry number matrix model of a biochemical network.However, most biochemical systems operate under non-equilibrium conditions.In fact, the true dissociation constant for a reaction is an empirical measure and is the ratio between the exponent (stoichiometry coefficients) forms for each reactant in a perfectly reversible reaction.The advent of large data warrants analytical approaches that depart significantly from previous methods in an effort to optimize computational resources whilst delivering meaningful information about the network.The use of fractional derivatives has gained traction over the last several years in an effort to impute missing data, improve convergence times and model complex systems [62][63][64][65][66] .For example, fractional derivatives are being used to improve convergence in the presence of multiple local optima and model and/or analyse the complex behaviour of fluids, biomolecules and biochemical systems such as anomalous diffusion in the cytoplasm and across membranes and even telomeres in the nucleoplasm [62][63][64][65][66] .All these studies have reported significant improvements in predicting an outcome, goodness-of-fit studies and time to convergence.In general, a fractional derivative of α-order can be approximated with the gamma function using the Rie- mann-Liouville (RL), Caputo (CAP) or the Grunwald-Letnikov (GL) methods [67][68][69][70] .Consider the Caputo fractional derivative, In the absence of empirical data investigating a biochemical network is challenging.One strategy for interrogating a biochemical network at non-equilibrium or near steady-state conditions is approximating the dissociation constant for a reaction.This equation is fundamental to biochemical analysis and involves a Figure 2. Computational studies with "ReDirection" to demonstrate suitability and biochemical relevance of the probable dissociation constant in characterizing a biochemical network for aerobic glycolysis.The dissociation constant for a reaction is an empirically derived parameter and can reveal valuable information such as the rate and rate constant for a reaction.Theoretically-derived approximations are usually based on user-defined bounds and/or pre-selecting strictly positive null space spanning vectors.Although the computations are rendered more feasible with these assumptions, the biomedical relevance of the data generated with these numeric constraints is likely to be of limited relevance.The algorithm presented computes the probable dissociation constant from first principles and with biochemically relevant constraints.A biochemical network of aerobic glycolysis is modelled and defined in terms of a sparse matrix of stoichiometric numbers of its reactants/products (J = 17) and reactions I = 20 ."ReDirection" checks this user-defined matrix and computes the null space, utilizes mathematically sound (tests of convergence, descriptive statistics, linear maps) constraints to compute the probable dissociation constant and thereby assign a dominant direction to every reaction of the modelled biochemical network of aerobic glycolysis.I , Total number of reactions of modelled biochemical network; J , Total number of reactants of modelled biochemical network; S p , reaction-centric and user-defined stoichiometry number matrix of a biochemical network for aerobic glycolysis; m, c; mitochondrial and cytosolic components of a reactant/product; η i , probable dissociation constant for the ith-reaction of a biochemical network with forward (F) , reverse (B) , or equivalent (E) outcomes.
phenomenological approximation of the integer values for the stoichiometry numbers of the reactants with an expression that involves real number exponents for each reactant of a fully reversible reaction.In accordance with Reaction 1 this is written as, We can then proceed to compute the fractional gradient descent, post hoc, to a preliminary imputation step to estimate missing or latent data points 62,63,69,70 .We can combine a small learning rate with a constant (stochastic) or mini-batch update protocol concomitantly with several novel algorithmic approaches to accomplish these steps 62,63,69,70 , The imputed values computed earlier will represent an initial approximation to the original stoichiometry numbers matrix for a biochemical network and can be refined iteratively, The computed numerical values can be used to assess the contribution of each molecular species to a reaction vector and are easily mapped to the exponents for each reactant.This can then be directly included in the formula to estimate the dissociation constant for a reaction at a particular step.Since this is a strictly positive real number we can use this metric to check whether each reaction vector converges to the expected outcome, where, where, where, www.nature.com/scientificreports/

Biomedical relevance of the presented algorithm to the study of biochemical networks
The dissociation constant for a reaction is a fundamental parameter in biochemical analysis and is computed predominantly from empirical data.However, this information is sparse and unavailable for several biochemical reactions of interest.For those reactions that laboratory data is available, the absence of corroborating in vivo details will render these less likely to have clinical-or biomedical-relevance. Additionally, biochemical assays comprise single or at most coupled reactions which further limits the utility of this information since complex biochemical function is brought about by biochemical networks.Although several algorithms utilize the null space, these are directed towards comprehending possible stable states for a modelled biochemical network.
Parameter selection is also left to the discretion of the user further limiting the portability of these studies.Unlike these investigations the probable dissociation constants are computed for every reaction and directly from a null space-generated subspace of the stoichiometry number matrix for a biochemical network.The algorithm concentrates on using combinatorial summations of unique and non-trivial null space spanning vectors to estimate the probable dissociation constant for every reaction along with mathematically rigorous and statistically relevant term-selection.A significant contribution of this algorithm is that it does not exclude numerical values and in doing so ensures that the computations are all encompassing.The computation of the probable dissociation constants are not only user-independent and computed from first principles, but may also be biochemically relevant.Since the mapping of the probable dissociation constants to the positive real number space is in accordance with established biochemical paradigms, the information gleaned can easily be incorporated into simulation studies for a biochemical network and used to compute the trajectories for a reactant/product 38,40 .Furthermore, the distribution of the probable dissociation constants across a biochemical network can offer plausible explanations into metabolite flux and usage.These testable hypotheses may be corroborated with detailed site-directed mutagenesis experiments and labelling studies among several others [44][45][46][47] .The probable dissociation constants can also be used to compare individual biochemical reactions of a network.These investigations may also be complemented by detailed theoretical studies such as phylogenetic analysis for enzyme-mediated reactions to ascertain the isoforms that may be deployed by a biochemical network.It may also be possible to perturb these networks and investigate the distribution for each pair of reactions of a biochemical network, i.e., before and after a stimulus.

Conclusions
The probable disassociation constant for a biochemical reaction is a numerical measure which can be computed directly from a null space-generated subspace of the stoichiometry number matrix for a biochemical network.Here, we present a mathematically rigorous algorithm to compute the probable dissociation constant for a reaction and thence for every reaction of a biochemical network.Each step of the algorithm is outlined and supported by the necessary mathematical formalism that is needed to comprehend and utilize this framework.The theoretical assertions presented in the main text are supported by formal proofs wherever applicable and are available and accessible as Supplementary Material.The biological relevance of this work is illustrated with computational studies and analyses of the probable dissociation constants for a well characterized biochemical network.The data generated is in accordance with established kinetic paradigms and is therefore, a suitable index of biochemical function and can be used to parameterize and study a biochemical network.These data are also amenable to detailed analyses and can be used to test various hypotheses for both, baseline-and perturbed-conditions of different cell types and across taxa.Additionally, the resulting data may offer valuable insights into the physiological and thence the pathological basis of diseases such as malignancies, altered immune response, etc.